/*
 * Copyright 2010-2012 Susanta Tewari. <freecode4susant@users.sourceforge.net>
 *
 * This program is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program.  If not, see <http://www.gnu.org/licenses/>.
 */

package coalescent.phylogeny;

import commons.util.ArrayUtil;
import commons.util.BinaryDataUtil;
import commons.util.ExceptionUtil;
import commons.util.IterableUtil;
import commons.util.graph.Edge;
import commons.util.graph.Node;

import java.util.*;

/**
 * An algorithm [1] by Dan Gusfield to test and build phylogeny for infinite-sites data. [1]
 * "Efficient Algorithms for Inferring Evolutionary Trees", By Dan Gusfield, Networks Vol. 21 (1991)
 * 19-28.
 *
 * @author Susanta Tewari
 * @since 1.3.0
 */
final class GusfieldAlgorithm {

    /** Field description */
    private final int[][] mPrime;


    // characters for mprime

    /** Field description */
    private final String[] mPrime_characters;

    /** Field description */
    private final int[][] LArray;

    /** Field description */
    private final int[] L;

    /** Field description */
    private final List<String> alleles;

    /** Field description */
    public final String NO_MUTATION_SYMBOL;

    /** Field description */
    public final String MUTATION_DELIM;

    /**
     * Computes intermediate structures of Gusfield`s algorithm to help in testing and building a
     * phylogeny from infinite-sites data.
     *
     * @param data infinite sites data (with segregating sites only ?); o indicates wild type and 1
     * indicates mutation
     * @param NO_MUTATION_SYMBOL a symbol indicating no mutation on an edge
     * @param MUTATION_DELIM delimeter to separate multiple mutations on an edge
     * @throws NullPointerException if any of the parameters is {@code null}
     * @throws IllegalArgumentException if, *) all arrays in {@code data} are not of same length *)
     * the arrays contain values other than 0 or 1
     */
    GusfieldAlgorithm(Iterable<byte[]> data, final String NO_MUTATION_SYMBOL,
                      final String MUTATION_DELIM) {

        if (NO_MUTATION_SYMBOL == null) {
            ExceptionUtil.throwIllegalArgExNull("NO_MUTATION_SYMBOL");
        }

        if (MUTATION_DELIM == null) {
            ExceptionUtil.throwIllegalArgExNull("MUTATION_DELIM");
        }


        // handles exceptions
        CheckPPhylogenyAlgo.checkData(data);

        this.NO_MUTATION_SYMBOL = NO_MUTATION_SYMBOL;
        this.MUTATION_DELIM     = MUTATION_DELIM;


        // create alleles (used for reporting errors)
        this.alleles = new ArrayList<>(10);

        for (byte[] val : data) {

            alleles.add(BinaryDataUtil.buildBinaryString(val));
        }


        // iterable to array
        byte[][] M            = IterableUtil.toArray(data);
        String[] m_characters = alleles.toArray(new String[0]);


        // STEP 1: Treat columns as binary numbers, considering the first row containing the most
        // significant bit. Sort these numbers in decresing order to sort the columns accordingly.
        // mPrime_characters would contain the sorted characters.
        // key: binary number value: column index, used in sorting column of M
        Map<Integer, Integer> indexMap    = new HashMap<>();
        Map<Integer, String> characterMap = new HashMap<>();
        String binaryString               = "";

        for (int j = 0; j < M[0].length; j++) {

            binaryString = "";

            for (byte[] aM : M) {

                binaryString += "" + aM[j];
            }

            Integer number = Integer.valueOf(binaryString, 2);

            if (characterMap.containsKey(number)) {

                characterMap.put(number,
                                 characterMap.get(number) + MUTATION_DELIM + m_characters[j]);

            } else {
                characterMap.put(number, "" + m_characters[j]);
            }

            indexMap.put(number, j);
        }


        // sorted map with decreasing order
        SortedMap<Integer, Integer> sortedMap = new TreeMap<>(Collections.reverseOrder());

        sortedMap.putAll(indexMap);

        Map<Integer, String> sortedCharMap = new TreeMap<>(Collections.reverseOrder());

        sortedCharMap.putAll(characterMap);

        this.mPrime_characters = sortedCharMap.values().toArray(new String[0]);


        // STEP 2: sort the columns of M, removing the duplicate columns, Call the new matrix M'
        mPrime = new int[M.length][sortedMap.size()];

        int counter = 0;

        for (Integer key : sortedMap.keySet()) {

            int colIndex = sortedMap.get(key);

            for (int i = 0; i < mPrime.length; i++) {

                mPrime[i][counter] = M[i][colIndex];
            }

            counter++;
        }


        // STEP 3: L[i][j] = Max (k) | k < j AND M'[i][k-1] == 1
        // = 0 if no such k exists
        // Here k is position rather than array index
        LArray = new int[mPrime.length][mPrime[0].length];
        L      = new int[mPrime[0].length];

        for (int j = 0; j < mPrime[0].length; j++) {

            int maxJ = 0;

            for (int i = 0; i < mPrime.length; i++) {

                int maxK = 0;

                for (int k = 0; k < j; k++) {

                    if (mPrime[i][k] == 1) {
                        maxK = k + 1;
                    }
                }

                LArray[i][j] = maxK;

                if ((mPrime[i][j] == 1) && (LArray[i][j] > maxJ)) {
                    maxJ = LArray[i][j];
                }
            }

            L[j] = maxJ;
        }
    }

    /**
     * Checks phylogeny of its data.
     *
     * @return {@code true} if it has phylogeny; {@code false} otherwise
     */
    public boolean isPhylogeny() {

        for (int i = 0; i < mPrime.length; i++) {

            for (int j = 0; j < mPrime[0].length; j++) {

                if ((mPrime[i][j] == 1) && (LArray[i][j] != L[j])) {
                    return false;
                }
            }
        }

        return true;
    }

    /**
     * @return the gene tree root or {@code null} if the data is not phylogeny
     */
    public Node<String, String> buildTree() {

        if (!isPhylogeny()) {
            return null;
        }


        // BUILD THE TREE
        // use this where the capacity is unknown or not efficient to compute
        final int UNKNOWN_CHILDREN_CAPACITY = 5;


        // create the root
        Node<String, String> root = new Node<>(null, "r", 10);


        // create the nodes and edges based on L()
        for (int i = 0; i < L.length; i++) {

            if (L[i] == 0) {
                root.addChild("n" + (i + 1), mPrime_characters[i], UNKNOWN_CHILDREN_CAPACITY);
            } else {

                Node<String, String> lookupNode = root.findFirstNode("n" + L[i]);

                lookupNode.addChild("n" + (i + 1), mPrime_characters[i], UNKNOWN_CHILDREN_CAPACITY);
            }
        }


        // create vector c using M'
        int[] c = new int[mPrime.length];

        for (int i = 0; i < mPrime.length; i++) {

            int position = 0;

            for (int j = 0; j < mPrime[0].length; j++) {

                if (mPrime[i][j] == 1) {
                    position = (j + 1);
                }
            }

            c[i] = position;
        }


        // process edges in c
        for (int i = 0; i < c.length; i++) {

            String nodeLabel = alleles.get(i);

            if (c[i] == 0) {

                root.addChild(nodeLabel, NO_MUTATION_SYMBOL, UNKNOWN_CHILDREN_CAPACITY);

                continue;
            }

            Edge<String, String> lookedupEdge = root.findFirstEdge(mPrime_characters[c[i] - 1]);
            Node<String, String> edgeHead     = lookedupEdge.getChild();

            if (edgeHead.isLeaf()) {
                edgeHead.setData(nodeLabel);
            } else {
                edgeHead.addChild(nodeLabel, NO_MUTATION_SYMBOL, UNKNOWN_CHILDREN_CAPACITY);
            }
        }

        return root;
    }

    /**
     * Method description
     *
     *
     * @return
     */
    public int[][] getMPrime() {
        return ArrayUtil.clone(mPrime);
    }

    /**
     * Method description
     *
     *
     * @return
     */
    public int[] getL() {
        return L.clone();
    }

    /**
     * Mutation labels sorted after the columns of M'. Mutation labels for the same columns are
     * concatenated {@link GusfieldAlgorithm#MUTATION_DELIM}.
     *
     * @return mutation characters (characters with the same data pattern are grouped)
     */
    public String[] getMPrimeCharacters() {
        return mPrime_characters.clone();
    }
}
